你的单细胞亚群可能并不是真的既有T细胞也有B细胞的标记基因

我们在前面的教程：为什么一个单细胞亚群既有T细胞也有B细胞的标记基因呢，给大家演示了一种可能性，就是cycling 亚群就是 既有T细胞也有B细胞的标记基因，实际上是因为它们这群细胞高表达 TOP2A和MKI67基因，处于比较活跃的细胞增殖状态，这样无论它本身是 T细胞还是B细胞，它们的这个状态使得它们在普通的降维聚类分群的流程里面就会被混合在一起。

在今年（2021）年初的文章《A single-cell map of intratumoral changes during anti-PD1 treatment of patients with breast cancer》，就提到了这样的现象，研究者们的解决方案是单独拿出来这个 cycling 亚群 进行继续细分，就可以得到原始的单细胞亚群了，因为在具体的细分亚群数据集里面，细胞周期影响全部的细胞所以它反而没有影响了。

实际上有一个更简单的检查方案，我们以《2021-GSE122930-心衰小鼠模型》为例子，虽然第7这个单细胞亚群里面的确实是可以看到既有T细胞也有B细胞的标记基因表达，如下所示：

既有T细胞也有B细胞的标记基因表达

但是这个第7群的离散度实在是太大了，甚至有可能仅仅是算法层面强行分群在一起了，就需要具体去看，它里面的全部的细胞是否同时表达T以及B细胞的基因，还是说， 仅仅是它们因为某种生物学特性（比如细胞周期）导致两个完全不同单细胞亚群被混合了。

但是这个第7群和第17群不一样，虽然17也是既有T细胞也有B细胞的标记基因表达，但是它的MKI67和TOP2A是高表达的，所以是因为细胞周期的影响。

我们首先看第7群

很简单的取子集的代码：

# 前面的 sce.all 是你自己对 GSE122930数据集进行降维聚类分群的seurat 对象 
sce=sce.all[,sce.all$CCA_snn_res.0.8==7]
sce
sce=CreateSeuratObject(sce@assays$RNA@counts)
sce <- NormalizeData(sce)
sce = FindVariableFeatures(sce)
sce = ScaleData(sce, vars.to.regress = c("nFeature_RNA", "percent_mito"))
sce = RunPCA(sce, npcs = 20)
sce = RunTSNE(sce, npcs = 20)
sce = RunUMAP(sce, dims = 1:10)
sce <- FindNeighbors(sce, dims = 1:10, verbose = FALSE)
sce <- FindClusters(sce, resolution = 0.5, verbose = FALSE)
DimPlot(sce) 

library(patchwork)
DoHeatmap(sce,features = cg)  + DimPlot(sce,label = T) 

library(stringr)
cg = str_to_title(
  c(
    'CD3D', 'CD3E', 'CD4','CD8A',
    'CD19', 'CD79A', 'MS4A1' 
  )
)
FeaturePlot(sce,features = cg)
DoHeatmap(sce,features = cg)
gl = list(
  Tcells =  cg[1:4],
  Bcells =  cg[5:7]
)
gl
AddModuleScore(object = sce,features = gl[1])
T_mat = sce@assays$RNA@counts[gl[[1]],]
B_mat = sce@assays$RNA@counts[gl[[2]],]

table(colSums(T_mat) > 5 ,
      colSums(B_mat) > 5 )

肉眼可以看到3基本上就是一个混合体，如图所示：

非常多的细胞居然真的是既有T细胞也有B细胞的标记基因表达

而且确实非常多的细胞居然真的是既有T细胞也有B细胞的标记基因表达，而且是同时表达啊：

> table(colSums(T_mat) > 5 ,
+       colSums(B_mat) > 5 )
       
        FALSE TRUE
  FALSE   210  372
  TRUE     64  122

但是，如果我们使用同样的代码，提取第17群，可以看到：

并不会 同时表达在同一个细胞里面

这个17群虽然是既有T细胞也有B细胞的标记基因表达，但它们其实并不会 同时表达在同一个细胞里面，仅仅是因为处于细胞增殖状态而被算法强行划分为了一个单细胞亚群而已！

那这个第7群到底是该如何解释呢？

什么样的情况，会导致这样的单细胞细胞亚群的出现，该文献：《Single Cell Sequencing of Mouse Heart Immune Infiltrate in Pressure Overload-Driven Heart Failure Reveals Extent of Immune Activation》，我看了看，原文似乎是没有这样的单细胞亚群出现：

 

如果我们以pbmc3k为例

同样的代码：

library(SeuratData) #加载seurat数据集  
getOption('timeout')
options(timeout=10000)
#InstallData("pbmc3k")  
data("pbmc3k")  
sce <- pbmc3k.final  
library(Seurat)

library(stringr)
cg = str_to_upper(
  c(
    'CD3D', 'CD3E', 'CD4','CD8A',
    'CD19', 'CD79A', 'MS4A1' 
  )
)
cg
FeaturePlot(sce,features = cg)
DoHeatmap(sce,features = cg,size = 3)
gl = list(
  Tcells =  cg[1:4],
  Bcells =  cg[5:7]
)
gl
AddModuleScore(object = sce,features = gl[1])
T_mat = sce@assays$RNA@counts[gl[[1]],]
B_mat = sce@assays$RNA@counts[gl[[2]],]

table(colSums(T_mat) > 1 ,
      colSums(B_mat) > 1 )

 

哪怕是如此严苛的条件：

> table(colSums(T_mat) > 1 ,
+       colSums(B_mat) > 1 )
       
        FALSE TRUE
  FALSE   913  321
  TRUE   1380   24

也仅仅是有24个细胞有疑似T和B细胞标记共表达。